Degradation statistics

baRulho:quantifying habitat-induced degradation of (animal) acoustic signals

Author

Marcelo Araya-Salas

Published

July 26, 2023

Source code, data and annotation protocol found at https://github.com/maRce10/barulho_paper

Purposes

  • Summarize degradation metrics on re-recorded signals from playback experiment at Bosque de Tlalpan, Mexico City, 2019

  • Run statistical analyses

 

Load packages

Code
source("https://raw.githubusercontent.com/maRce10/sketchy/main/R/load_packages.R")

# install/ load packages
load_packages(packages = c("remotes", "kableExtra", "knitr", "formatR",
    "rprojroot", "maRce10/warbleR", "ggplot2", "tidyr", "viridis",
    "corrplot", "brms", "ggdist", "cowplot", "cmdstanr", "maRce10/brmsish",
    "emmeans", "ggsignif"))

my.viridis <- function(...) viridis(alpha = 0.5, begin = 0.3, end = 0.7,
    ...)

source("~/Dropbox/R_package_testing/brmsish/R/extended_summary.R")
source("~/Dropbox/R_package_testing/brmsish/R/helpers.R")
source("~/Dropbox/R_package_testing/brmsish/R/check_rds_fits.R")
Code
degrad_df <- read.csv("./data/processed/tlalpan_degradation_metrics_v01.csv")

degrad_params <- c("blur.ratio", "spectrum.blur.ratio", "envelope.correlation",
    "excess.attenuation", "signal.to.noise.ratio", "cross.correlation",
    "tail.to.signal.ratio", "tail.to.noise.ratio", "spectrum.correlation")

comp.cases <- complete.cases(degrad_df[, names(degrad_df) %in% degrad_params])

pca <- prcomp(degrad_df[comp.cases, names(degrad_df) %in% degrad_params],
    scale. = TRUE)

# add to data
degrad_df$pc1 <- NA
degrad_df$pc1[comp.cases] <- pca$x[, 1]
degrad_df$pc1.1m.rate <- degrad_df$distance
# plot rotation values by PC
pca_rot <- as.data.frame(pca$rotation[, 1:4])

pca_rot_stck <- stack(pca_rot)

pca_rot_stck$variable <- rownames(pca_rot)
pca_rot_stck$Sign <- ifelse(pca_rot_stck$values > 0, "Positive", "Negative")
pca_rot_stck$rotation <- abs(pca_rot_stck$values)

ggplot(pca_rot_stck, aes(x = variable, y = rotation, fill = Sign)) +
    geom_col() + coord_flip() + scale_fill_viridis_d(alpha = 0.7,
    begin = 0.2, end = 0.8) + facet_wrap(~ind) + theme_classic()

1 Correlation among metrics

Raw metrics:

Code
cormat <- cor(degrad_df[, degrad_params], use = "pairwise.complete.obs")

rownames(cormat) <- colnames(cormat) <- degrad_params

cols_corr <- colorRampPalette(c(viridis(3, direction = 1, begin = 0.2,
    end = 0.5), "#BEBEBE1A", "white", "#BEBEBE1A", viridis(3, direction = 1,
    begin = 0.7, end = 0.9)))(30)

cp <- corrplot.mixed(cormat, tl.cex = 0.7, upper.col = cols_corr,
    lower.col = cols_corr, order = "hclust", lower = "number", upper = "ellipse",
    tl.col = "black")

Code
# sort parameters as in clusters for cross correlation
degrad_params <- degrad_params[match(rownames(cp$corr), degrad_params)]

2 Data description

  • 20 frequencies
  • 3 locations
  • 11520 test sounds
  • 160 treatment combinations

Sample sizes per location, transect and signal type:

Code
agg <- aggregate(cbind(sound.id, treatment.replicates) ~ location +
    habitat.structure + distance, degrad_df, function(x) length(unique(x)))

agg$replicates <- agg$sound.id/agg$treatment.replicates

df1 <- knitr::kable(agg, row.names = FALSE, escape = FALSE, format = "html")

kable_styling(df1, bootstrap_options = c("striped", "hover", "condensed",
    "responsive"), full_width = FALSE, font_size = 15)
location habitat.structure distance sound.id treatment.replicates replicates
1 closed 10 480 160 3
2 closed 10 480 160 3
3 closed 10 480 160 3
1 open 10 480 160 3
2 open 10 480 160 3
3 open 10 480 160 3
1 closed 30 480 160 3
2 closed 30 480 160 3
3 closed 30 480 160 3
1 open 30 480 160 3
2 open 30 480 160 3
3 open 30 480 160 3
1 closed 65 480 160 3
2 closed 65 480 160 3
3 closed 65 480 160 3
1 open 65 480 160 3
2 open 65 480 160 3
3 open 65 480 160 3
1 closed 100 480 160 3
2 closed 100 480 160 3
3 closed 100 480 160 3
1 open 100 480 160 3
2 open 100 480 160 3
3 open 100 480 160 3

3 Statistical analysis

Code
iter <- 20000
chains <- 4
priors <- c(prior(normal(0, 6), class = "b"), prior(cauchy(0, 4),
    class = "sd"))

# set frequency to mean-centered and divide by twice the
# standard deviation
degrad_df$raw_frequency <- degrad_df$frequency
degrad_df$frequency <- (degrad_df$frequency - mean(degrad_df$frequency))/(2 *
    sd(degrad_df$frequency))

# set base level for factors
degrad_df$habitat.structure <- factor(degrad_df$habitat.structure,
    levels = c("open", "closed"))
degrad_df$frequency.modulation <- factor(degrad_df$frequency.modulation,
    levels = c("no_fm", "fm"))
degrad_df$amplitude.modulation <- factor(degrad_df$amplitude.modulation,
    levels = c("no_am", "am"))
degrad_df$duration <- factor(degrad_df$duration, levels = c("short",
    "long"))
degrad_df$location <- as.factor(degrad_df$location)

degrad_df$distance_f <- paste0(degrad_df$distance, "m")
degrad_df$distance_f <- factor(degrad_df$distance_f, levels = c("10m",
    "30m", "65m", "100m"), ordered = TRUE)


set.seed(123)

cmdstanr::set_cmdstan_path("~/Documentos/cmdstan/")

# to run within-chain parallelization
mod_pc1 <- brm(formula = pc1 ~ frequency * habitat.structure + frequency.modulation *
    habitat.structure + amplitude.modulation * habitat.structure +
    duration * habitat.structure + mo(distance_f) + (1 | location) +
    (1 | treatment.replicates), data = degrad_df, prior = priors,
    iter = iter, chains = chains, cores = chains, control = list(adapt_delta = 0.99,
        max_treedepth = 15), file = "./data/processed/regression_model_pc1.RDS",
    file_refit = "always", backend = "cmdstanr", threads = threading(10))

mod_blurratio <- brm(formula = blur.ratio ~ frequency * habitat.structure +
    frequency.modulation * habitat.structure + amplitude.modulation *
    habitat.structure + duration * habitat.structure + mo(distance_f) +
    (1 | location) + (1 | treatment.replicates), data = degrad_df,
    prior = priors, iter = iter, chains = chains, cores = chains,
    control = list(adapt_delta = 0.99, max_treedepth = 15), file = "./data/processed/regression_model_blur_ratio.RDS",
    file_refit = "always", backend = "cmdstanr", threads = threading(10))

mod_ea <- brm(formula = excess.attenuation ~ frequency * habitat.structure +
    frequency.modulation * habitat.structure + amplitude.modulation *
    habitat.structure + duration * habitat.structure + mo(distance_f) +
    (1 | location) + (1 | treatment.replicates), data = degrad_df,
    prior = priors, iter = iter, chains = chains, cores = chains,
    control = list(adapt_delta = 0.99, max_treedepth = 15), file = "./data/processed/regression_model_excess_attenuation.RDS",
    file_refit = "always", backend = "cmdstanr", threads = threading(10))

mod_tsr <- brm(formula = tail.to.signal.ratio ~ frequency * habitat.structure +
    frequency.modulation * habitat.structure + amplitude.modulation *
    habitat.structure + duration * habitat.structure + mo(distance_f) +
    (1 | location) + (1 | treatment.replicates), data = degrad_df,
    prior = priors, iter = iter, chains = chains, cores = chains,
    control = list(adapt_delta = 0.99, max_treedepth = 15), file = "./data/processed/regression_model_tail_to_signal_ratio.RDS",
    file_refit = "always", backend = "cmdstanr", threads = threading(10))

3.1 PC1 degradation

Code
effects <- c("habitat structure", "frequency modulation", "amplitude modulation",
    "frequency", "duration", "frequency:habitat structure", "habitat structure:duration",
    "habitat structure:amplitude modulation", "habitat structure:frequency modulation")

mod <- readRDS("./data/processed/regression_model_pc1.RDS")

extended_summary(mod, n.posterior = 1000, fill = viridis(10)[7], trace.palette = my.viridis,
    remove.intercepts = TRUE, highlight = TRUE, gsub.pattern = c("b_",
        "habitat.structureclosed", "amplitude.modulationam", "frequency.modulationfm",
        "durationlong"), gsub.replacement = c("", "habitat structure",
        "amplitude modulation", "frequency modulation", "duration"),
    effects = effects, print.name = FALSE, trace = TRUE, return = FALSE)
priors formula iterations chains thinning warmup diverg_transitions rhats > 1.05 min_bulk_ESS min_tail_ESS seed
1 b-normal(0, 6) Intercept-student_t(3, -0.3, 2.6) sd-cauchy(0, 4) sigma-student_t(3, 0, 2.6) simo-dirichlet(1) pc1 ~ frequency * habitat.structure + frequency.modulation * habitat.structure + amplitude.modulation * habitat.structure + duration * habitat.structure + mo(distance_f) + (1 | location) + (1 | treatment.replicates) 20000 4 1 10000 166 (0.00415%) 0 5197.53 5998.95 721775339
Estimate l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
habitat structure 1.064 0.981 1.148 1 26827.620 27957.994
frequency modulation 0.810 0.647 0.978 1 5503.606 10258.486
amplitude modulation 0.234 0.066 0.399 1.001 5271.999 9438.670
frequency 0.894 0.725 1.058 1 5274.519 9535.988
duration 0.151 -0.013 0.317 1.001 5197.530 10390.299
frequency:habitat structure 0.320 0.236 0.403 1 40733.644 26946.792
habitat structure:duration 0.081 -0.003 0.164 1 35912.172 27252.667
habitat structure:amplitude modulation -0.018 -0.101 0.066 1 28829.607 5998.950
habitat structure:frequency modulation -0.154 -0.237 -0.070 1 37758.911 29477.077

3.1.1 Contrasts

3.1.1.1 Frequency modulation * habitat structure

Code
cntrs <- as.data.frame(emmeans(mod, pairwise ~ frequency.modulation *
    habitat.structure)$contrasts)
names(cntrs)[3:4] <- c("l-95% CI", "u-95% CI")

coef_table <- html_format_coef_table(cntrs, highlight = TRUE, fill = viridis(10)[7])

print(coef_table)
contrast estimate l-95% CI u-95% CI
1 no_fm open - fm open -0.810 -0.977 -0.646
2 no_fm open - no_fm closed -1.096 -1.152 -1.034
3 no_fm open - fm closed -1.752 -1.915 -1.584
4 fm open - no_fm closed -0.285 -0.452 -0.122
5 fm open - fm closed -0.942 -1.000 -0.884
6 no_fm closed - fm closed -0.657 -0.821 -0.489

3.1.1.2 Amplitude modulation * habitat structure

Code
cntrs <- as.data.frame(emmeans(mod, pairwise ~ amplitude.modulation *
    habitat.structure)$contrasts)
names(cntrs)[3:4] <- c("l-95% CI", "u-95% CI")

coef_table <- html_format_coef_table(cntrs, highlight = TRUE, fill = viridis(10)[7])

print(coef_table)
contrast estimate l-95% CI u-95% CI
1 no_am open - am open -0.234 -0.394 -0.061
2 no_am open - no_am closed -1.028 -1.087 -0.970
3 no_am open - am closed -1.243 -1.407 -1.071
4 am open - no_am closed -0.794 -0.958 -0.623
5 am open - am closed -1.010 -1.069 -0.952
6 no_am closed - am closed -0.216 -0.383 -0.048

3.1.2 Conditional plots

Code
bs <- 10
ts <- 4
ylm <- c(-5, 1.5)
diff <- 0.1
sd_freq <- sd(degrad_df$frequency)
mean_freq <- round(mean(degrad_df$frequency))

freq_labs <- c(paste(round(mean_freq + sd_freq, 1), "(mean+SD)", sep = "\n"),
    paste(round(mean_freq - sd_freq, 1), "(mean-SD)", sep = "\n"))

ggf <- plot(conditional_effects(mod, "habitat.structure:frequency",
    int_conditions = list(frequency = c(mean_freq - sd_freq, mean_freq +
        sd_freq))), plot = FALSE)[[1]] + scale_color_viridis_d(end = 0.8,
    labels = freq_labs) + theme_classic(base_size = bs) + labs(color = "Frequency",
    y = "Degradation (PC1)", x = "Habitat") + scale_fill_discrete(guide = "none") +
    ylim(ylm)


# add contrasts
for (i in 1:2) ggf <- ggf + geom_signif(y_position = c(1), xmin = c(i -
    diff), xmax = c(i + diff), annotation = c("*"), tip_length = 0.02,
    textsize = ts, size = 0.4, color = "gray40")

# add larger bracket
ggf <- ggf + geom_signif(y_position = c(1.5), xmin = c(1), xmax = c(2),
    annotation = c("*"), tip_length = 0.02, textsize = ts, size = 0.4,
    color = "gray40")


ggfm <- plot(conditional_effects(mod, "habitat.structure:frequency.modulation"),
    plot = FALSE)[[1]] + scale_color_viridis_d(end = 0.8, labels = c("Tonal",
    "FM")) + theme_classic(base_size = bs) + labs(x = "Habitat", y = "",
    color = "") + scale_fill_discrete(guide = "none") + theme(legend.position = c(0.93,
    0.5)) + ylim(ylm)


# add contrasts
for (i in 1:2) ggfm <- ggfm + geom_signif(y_position = c(0), xmin = c(i -
    diff), xmax = c(i + diff), annotation = c("*"), tip_length = 0.02,
    textsize = ts, size = 0.4, color = "gray40")

ggfm <- ggfm + geom_signif(y_position = c(1), xmin = 1, xmax = 2,
    annotation = c("*"), tip_length = 0.02, textsize = ts, size = 0.4,
    color = "gray40")

ggam <- plot(conditional_effects(mod, "habitat.structure:amplitude.modulation"),
    plot = FALSE)[[1]] + scale_color_viridis_d(end = 0.8, labels = c("Flat",
    "AM")) + theme_classic(base_size = bs) + labs(x = "Habitat", y = "Degradation (PC1)",
    color = "") + scale_fill_discrete(guide = "none") + theme(legend.position = c(0.93,
    0.5)) + ylim(ylm)

# add contrasts
for (i in 1:2) ggam <- ggam + geom_signif(y_position = c(0), xmin = c(i -
    diff), xmax = c(i + diff), annotation = c("*"), tip_length = 0.02,
    textsize = ts, size = 0.4, color = "gray40")

ggam <- ggam + geom_signif(y_position = c(1), xmin = 1, xmax = 2,
    annotation = c("ns"), tip_length = 0.02, textsize = ts, size = 0.4,
    color = "gray40")


ggdur <- plot(conditional_effects(mod, "habitat.structure:duration"),
    plot = FALSE)[[1]] + scale_color_viridis_d(end = 0.8, labels = c("Short",
    "Long")) + theme_classic(base_size = bs) + labs(x = "Habitat",
    y = "", color = "") + scale_fill_discrete(guide = "none") + theme(legend.position = c(0.93,
    0.5)) + ylim(ylm)

# add contrasts
for (i in 1:2) ggdur <- ggdur + geom_signif(y_position = c(0), xmin = c(i -
    diff), xmax = c(i + diff), annotation = c("ns"), tip_length = 0.02,
    textsize = ts, size = 0.4, color = "gray40")

ggdur <- ggdur + geom_signif(y_position = c(1), xmin = 1, xmax = 2,
    annotation = c("ns"), tip_length = 0.02, textsize = 4, size = 0.4,
    color = "gray40")


plot_grid(ggf, ggfm, ggam, ggdur)

3.2 Blur ratio

Code
mod <- readRDS("./data/processed/regression_model_blur_ratio.RDS")

extended_summary(fit = mod, n.posterior = 1000, fill = viridis(10)[7],
    trace.palette = my.viridis, remove.intercepts = TRUE, highlight = TRUE,
    gsub.pattern = c("b_", "habitat.structureclosed", "amplitude.modulationam",
        "frequency.modulationfm", "durationlong"), gsub.replacement = c("",
        "habitat structure", "amplitude modulation", "frequency modulation",
        "duration"), effects = effects, print.name = FALSE)
priors formula iterations chains thinning warmup diverg_transitions rhats > 1.05 min_bulk_ESS min_tail_ESS seed
1 b-normal(0, 6) Intercept-student_t(3, 0.1, 2.5) sd-cauchy(0, 4) sigma-student_t(3, 0, 2.5) simo-dirichlet(1) blur.ratio ~ frequency * habitat.structure + frequency.modulation * habitat.structure + amplitude.modulation * habitat.structure + duration * habitat.structure + mo(distance_f) + (1 | location) + (1 | treatment.replicates) 20000 4 1 10000 220 (0.0055%) 0 3942.862 7547.377 1717938324
Estimate l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
habitat structure 0.029 0.026 0.033 1 28754.286 28266.622
frequency modulation 0.065 0.058 0.072 1 4252.830 8903.955
amplitude modulation 0.024 0.017 0.031 1 4657.769 9172.275
frequency 0.007 0.000 0.014 1.001 4108.394 7547.377
duration 0.001 -0.006 0.008 1.001 3942.862 8754.354
frequency:habitat structure 0.016 0.013 0.020 1 52057.556 28431.402
habitat structure:duration 0.005 0.001 0.009 1 41292.428 31201.185
habitat structure:amplitude modulation 0.002 -0.001 0.006 1 40019.432 30015.581
habitat structure:frequency modulation -0.020 -0.024 -0.017 1 40647.068 30954.625

3.2.1 Contrasts

3.2.1.1 Frequency modulation * habitat structure

Code
cntrs <- as.data.frame(emmeans(mod, pairwise ~ frequency.modulation *
    habitat.structure)$contrasts)
names(cntrs)[3:4] <- c("l-95% CI", "u-95% CI")

coef_table <- html_format_coef_table(cntrs, highlight = TRUE, fill = viridis(10)[7])

print(coef_table)
contrast estimate l-95% CI u-95% CI
1 no_fm open - fm open -0.065 -0.072 -0.058
2 no_fm open - no_fm closed -0.033 -0.035 -0.030
3 no_fm open - fm closed -0.078 -0.084 -0.070
4 fm open - no_fm closed 0.033 0.026 0.040
5 fm open - fm closed -0.012 -0.015 -0.010
6 no_fm closed - fm closed -0.045 -0.052 -0.038

3.2.1.2 Amplitude modulation * habitat structure

Code
cntrs <- as.data.frame(emmeans(mod, pairwise ~ amplitude.modulation *
    habitat.structure)$contrasts)
names(cntrs)[3:4] <- c("l-95% CI", "u-95% CI")

coef_table <- html_format_coef_table(cntrs, highlight = TRUE, fill = viridis(10)[7])

print(coef_table)
contrast estimate l-95% CI u-95% CI
1 no_am open - am open -0.024 -0.031 -0.017
2 no_am open - no_am closed -0.021 -0.024 -0.019
3 no_am open - am closed -0.047 -0.054 -0.040
4 am open - no_am closed 0.002 -0.005 0.009
5 am open - am closed -0.024 -0.026 -0.021
6 no_am closed - am closed -0.026 -0.033 -0.019

3.2.2 Conditional plots

Code
bs <- 10
ts <- 4
ylim <- c(-0.1, 0.3)

ggf <- plot(conditional_effects(mod, "habitat.structure:frequency",
    int_conditions = list(frequency = c(mean_freq - sd_freq, mean_freq +
        sd_freq))), plot = FALSE)[[1]] + scale_color_viridis_d(end = 0.8,
    labels = freq_labs) + theme_classic(base_size = bs) + labs(color = "Frequency",
    y = "Blur ratio", x = "Habitat") + scale_fill_discrete(guide = "none") +
    ylim(ylim)

# add contrasts
for (i in 1:2) ggf <- ggf + geom_signif(y_position = c(0.2), xmin = c(i -
    diff), xmax = c(i + diff), annotation = c("*"), tip_length = 0.02,
    textsize = ts, size = 0.4, color = "gray40")

# add larger bracket
ggf <- ggf + geom_signif(y_position = c(0.23), xmin = c(1), xmax = c(2),
    annotation = c("*"), tip_length = 0.02, textsize = ts, size = 0.4,
    color = "gray40")



ggfm <- plot(conditional_effects(mod, "habitat.structure:frequency.modulation"),
    plot = FALSE)[[1]] + scale_color_viridis_d(end = 0.8, labels = c("Tonal",
    "FM")) + theme_classic(base_size = bs) + labs(x = "Habitat", y = "",
    color = "") + scale_fill_discrete(guide = "none") + theme(legend.position = c(0.93,
    0.5)) + ylim(ylim)

diff <- 0.1


# add contrasts
for (i in 1:2) ggfm <- ggfm + geom_signif(y_position = c(0.25), xmin = c(i -
    diff), xmax = c(i + diff), annotation = c("*"), tip_length = 0.02,
    textsize = ts, size = 0.4, color = "gray40")

ggfm <- ggfm + geom_signif(y_position = c(0.3), xmin = 1, xmax = 2,
    annotation = c("*"), tip_length = 0.02, textsize = ts, size = 0.4,
    color = "gray40")

ggam <- plot(conditional_effects(mod, "habitat.structure:amplitude.modulation"),
    plot = FALSE)[[1]] + scale_color_viridis_d(end = 0.8, labels = c("Flat",
    "AM")) + theme_classic(base_size = bs) + labs(x = "Habitat", y = "Blur ratio",
    color = "") + scale_fill_discrete(guide = "none") + theme(legend.position = c(0.93,
    0.5)) + ylim(ylim)

# add contrasts
for (i in 1:2) ggam <- ggam + geom_signif(y_position = c(0.2), xmin = c(i -
    diff), xmax = c(i + diff), annotation = c("*"), tip_length = 0.02,
    textsize = ts, size = 0.4, color = "gray40")

ggam <- ggam + geom_signif(y_position = c(0.3), xmin = 1, xmax = 2,
    annotation = c("ns"), tip_length = 0.02, textsize = ts, size = 0.4,
    color = "gray40")


ggdur <- plot(conditional_effects(mod, "habitat.structure:duration"),
    plot = FALSE)[[1]] + scale_color_viridis_d(end = 0.8, labels = c("Short",
    "Long")) + theme_classic(base_size = bs) + labs(x = "Habitat",
    y = "", color = "") + scale_fill_discrete(guide = "none") + theme(legend.position = c(0.93,
    0.5)) + ylim(ylim)

# add contrasts
for (i in 1:2) ggdur <- ggdur + geom_signif(y_position = c(0.2), xmin = c(i -
    diff), xmax = c(i + diff), annotation = c("ns"), tip_length = 0.02,
    textsize = ts, size = 0.4, color = "gray40")

ggdur <- ggdur + geom_signif(y_position = c(0.3), xmin = 1, xmax = 2,
    annotation = c("ns"), tip_length = 0.02, textsize = 4, size = 0.4,
    color = "gray40")


plot_grid(ggf, ggfm, ggam, ggdur)

3.3 Excess attenuation

Code
mod <- readRDS("./data/processed/regression_model_excess_attenuation.RDS")

extended_summary(fit = mod, n.posterior = 1000, fill = viridis(10)[7],
    trace.palette = my.viridis, remove.intercepts = TRUE, highlight = TRUE,
    gsub.pattern = c("b_", "habitat.structureclosed", "amplitude.modulationam",
        "frequency.modulationfm", "durationlong"), gsub.replacement = c("",
        "habitat structure", "amplitude modulation", "frequency modulation",
        "duration"), effects = effects, print.name = FALSE)
priors formula iterations chains thinning warmup diverg_transitions rhats > 1.05 min_bulk_ESS min_tail_ESS seed
1 b-normal(0, 6) Intercept-student_t(3, 83.8, 34) sd-cauchy(0, 4) sigma-student_t(3, 0, 34) simo-dirichlet(1) excess.attenuation ~ frequency * habitat.structure + frequency.modulation * habitat.structure + amplitude.modulation * habitat.structure + duration * habitat.structure + mo(distance_f) + (1 | location) + (1 | treatment.replicates) 20000 4 1 10000 18 (0.00045%) 0 4940.947 8799.022 1879993832
Estimate l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
habitat structure 9.891 9.055 10.729 1 24779.171 28856.313
frequency modulation -3.169 -4.685 -1.650 1 4995.106 10015.215
amplitude modulation -2.694 -4.235 -1.143 1.001 4978.640 9873.109
frequency 17.937 16.401 19.482 1 4940.947 8799.022
duration -0.848 -2.374 0.693 1.001 4990.889 9951.945
frequency:habitat structure -3.621 -4.464 -2.781 1 52750.860 29230.209
habitat structure:duration -0.399 -1.230 0.445 1 39494.566 31368.784
habitat structure:amplitude modulation -1.574 -2.411 -0.739 1 37511.727 30680.835
habitat structure:frequency modulation -0.131 -0.978 0.713 1 38959.256 31041.548

3.3.1 Contrasts

3.3.1.1 Frequency modulation * habitat structure

Code
cntrs <- as.data.frame(emmeans(mod, pairwise ~ frequency.modulation *
    habitat.structure)$contrasts)
names(cntrs)[3:4] <- c("l-95% CI", "u-95% CI")

coef_table <- html_format_coef_table(cntrs, highlight = TRUE, fill = viridis(10)[7])

print(coef_table)
contrast estimate l-95% CI u-95% CI
1 no_fm open - fm open 3.173 1.607 4.639
2 no_fm open - no_fm closed -8.903 -9.512 -8.319
3 no_fm open - fm closed -5.602 -7.118 -4.080
4 fm open - no_fm closed -12.073 -13.596 -10.569
5 fm open - fm closed -8.772 -9.376 -8.179
6 no_fm closed - fm closed 3.302 1.763 4.815

3.3.1.2 Amplitude modulation * habitat structure

Code
cntrs <- as.data.frame(emmeans(mod, pairwise ~ amplitude.modulation *
    habitat.structure)$contrasts)
names(cntrs)[3:4] <- c("l-95% CI", "u-95% CI")

coef_table <- html_format_coef_table(cntrs, highlight = TRUE, fill = viridis(10)[7])

print(coef_table)
contrast estimate l-95% CI u-95% CI
1 no_am open - am open 2.696 1.095 4.181
2 no_am open - no_am closed -9.625 -10.207 -9.023
3 no_am open - am closed -5.360 -6.850 -3.789
4 am open - no_am closed -12.325 -13.895 -10.811
5 am open - am closed -8.053 -8.647 -7.454
6 no_am closed - am closed 4.270 2.704 5.761

3.3.2 Conditional plots

Code
bs <- 10
ts <- 4
ylim <- c(20, 90)

ggf <- plot(conditional_effects(mod, "habitat.structure:frequency",
    int_conditions = list(frequency = c(mean_freq - sd_freq, mean_freq +
        sd_freq))), plot = FALSE)[[1]] + scale_color_viridis_d(end = 0.8,
    labels = freq_labs) + theme_classic(base_size = bs) + labs(color = "Frequency",
    y = "Excess attenuation", x = "Habitat") + scale_fill_discrete(guide = "none") +
    ylim(ylim)

# add contrasts
for (i in 1:2) ggf <- ggf + geom_signif(y_position = c(85), xmin = c(i -
    diff), xmax = c(i + diff), annotation = c("*"), tip_length = 0.02,
    textsize = ts, size = 0.4, color = "gray40")

# add larger bracket
ggf <- ggf + geom_signif(y_position = c(90), xmin = c(1), xmax = c(2),
    annotation = c("*"), tip_length = 0.02, textsize = ts, size = 0.4,
    color = "gray40")



ggfm <- plot(conditional_effects(mod, "habitat.structure:frequency.modulation"),
    plot = FALSE)[[1]] + scale_color_viridis_d(end = 0.8, labels = c("Tonal",
    "FM")) + theme_classic(base_size = bs) + labs(x = "Habitat", y = "",
    color = "") + scale_fill_discrete(guide = "none") + theme(legend.position = c(0.93,
    0.5)) + ylim(ylim)

diff <- 0.1


# add contrasts
for (i in 1:2) ggfm <- ggfm + geom_signif(y_position = c(65), xmin = c(i -
    diff), xmax = c(i + diff), annotation = c("*"), tip_length = 0.02,
    textsize = ts, size = 0.4, color = "gray40")

ggfm <- ggfm + geom_signif(y_position = c(72), xmin = 1, xmax = 2,
    annotation = c("ns"), tip_length = 0.02, textsize = ts, size = 0.4,
    color = "gray40")

ggam <- plot(conditional_effects(mod, "habitat.structure:amplitude.modulation"),
    plot = FALSE)[[1]] + scale_color_viridis_d(end = 0.8, labels = c("Flat",
    "AM")) + theme_classic(base_size = bs) + labs(x = "Habitat", y = "Excess attenuation",
    color = "") + scale_fill_discrete(guide = "none") + theme(legend.position = c(0.93,
    0.5)) + ylim(ylim)

# add contrasts
for (i in 1:2) ggam <- ggam + geom_signif(y_position = c(65), xmin = c(i -
    diff), xmax = c(i + diff), annotation = c("*"), tip_length = 0.02,
    textsize = ts, size = 0.4, color = "gray40")

ggam <- ggam + geom_signif(y_position = c(72), xmin = 1, xmax = 2,
    annotation = c("*"), tip_length = 0.02, textsize = ts, size = 0.4,
    color = "gray40")


ggdur <- plot(conditional_effects(mod, "habitat.structure:duration"),
    plot = FALSE)[[1]] + scale_color_viridis_d(end = 0.8, labels = c("Short",
    "Long")) + theme_classic(base_size = bs) + labs(x = "Habitat",
    y = "", color = "") + scale_fill_discrete(guide = "none") + theme(legend.position = c(0.93,
    0.5)) + ylim(ylim)

# add contrasts
for (i in 1:2) ggdur <- ggdur + geom_signif(y_position = c(65), xmin = c(i -
    diff), xmax = c(i + diff), annotation = c("ns"), tip_length = 0.02,
    textsize = ts, size = 0.4, color = "gray40")

ggdur <- ggdur + geom_signif(y_position = c(72), xmin = 1, xmax = 2,
    annotation = c("ns"), tip_length = 0.02, textsize = 4, size = 0.4,
    color = "gray40")


plot_grid(ggf, ggfm, ggam, ggdur)

3.4 Tail-to-signal ratio

Code
mod <- readRDS("./data/processed/regression_model_tail_to_signal_ratio.RDS")

extended_summary(fit = mod, n.posterior = 1000, fill = viridis(10)[7],
    trace.palette = my.viridis, remove.intercepts = TRUE, highlight = TRUE,
    gsub.pattern = c("b_", "habitat.structureclosed", "amplitude.modulationam",
        "frequency.modulationfm", "durationlong"), gsub.replacement = c("",
        "habitat structure", "amplitude modulation", "frequency modulation",
        "duration"), effects = effects, print.name = FALSE)
priors formula iterations chains thinning warmup diverg_transitions rhats > 1.05 min_bulk_ESS min_tail_ESS seed
1 b-normal(0, 6) Intercept-student_t(3, -6.3, 8.7) sd-cauchy(0, 4) sigma-student_t(3, 0, 8.7) simo-dirichlet(1) tail.to.signal.ratio ~ frequency * habitat.structure + frequency.modulation * habitat.structure + amplitude.modulation * habitat.structure + duration * habitat.structure + mo(distance_f) + (1 | location) + (1 | treatment.replicates) 20000 4 1 10000 18 (0.00045%) 0 4589.723 9302.75 1096594823
Estimate l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
habitat structure 3.254 2.899 3.613 1 22477.025 27029.479
frequency modulation -0.815 -1.436 -0.195 1.001 4589.723 9368.327
amplitude modulation 0.822 0.218 1.432 1.001 4978.120 10414.482
frequency 1.746 1.138 2.362 1.001 4828.315 9302.750
duration -0.426 -1.040 0.187 1.001 4863.510 9871.931
frequency:habitat structure 0.116 -0.235 0.469 1 43544.100 29786.404
habitat structure:duration -0.056 -0.410 0.296 1 34372.988 30274.700
habitat structure:amplitude modulation -0.384 -0.739 -0.029 1 33163.575 29602.429
habitat structure:frequency modulation 0.539 0.181 0.888 1 32846.121 31206.849

Code
# mod <- readRDS('./data/processed/regression_model_int.RDS')

3.4.1 Contrasts

3.4.1.1 Frequency modulation * habitat structure

Code
cntrs <- as.data.frame(emmeans(mod, pairwise ~ frequency.modulation *
    habitat.structure)$contrasts)
names(cntrs)[3:4] <- c("l-95% CI", "u-95% CI")

coef_table <- html_format_coef_table(cntrs, highlight = TRUE, fill = viridis(10)[7])

print(coef_table)
contrast estimate l-95% CI u-95% CI
1 no_fm open - fm open 0.815 0.205 1.443
2 no_fm open - no_fm closed -3.034 -3.284 -2.788
3 no_fm open - fm closed -2.759 -3.372 -2.133
4 fm open - no_fm closed -3.849 -4.463 -3.228
5 fm open - fm closed -3.573 -3.828 -3.324
6 no_fm closed - fm closed 0.275 -0.343 0.894

3.4.1.2 Amplitude modulation * habitat structure

Code
cntrs <- as.data.frame(emmeans(mod, pairwise ~ amplitude.modulation *
    habitat.structure)$contrasts)
names(cntrs)[3:4] <- c("l-95% CI", "u-95% CI")

coef_table <- html_format_coef_table(cntrs, highlight = TRUE, fill = viridis(10)[7])

print(coef_table)
contrast estimate l-95% CI u-95% CI
1 no_am open - am open -0.821 -1.444 -0.230
2 no_am open - no_am closed -3.495 -3.749 -3.247
3 no_am open - am closed -3.933 -4.546 -3.329
4 am open - no_am closed -2.676 -3.284 -2.065
5 am open - am closed -3.111 -3.366 -2.862
6 no_am closed - am closed -0.437 -1.062 0.153

3.4.2 Conditional plots

Code
bs <- 10
ts <- 4
ylim <- c(-27, -5)

ggf <- plot(conditional_effects(mod, "habitat.structure:frequency",
    int_conditions = list(frequency = c(mean_freq - sd_freq, mean_freq +
        sd_freq))), plot = FALSE)[[1]] + scale_color_viridis_d(end = 0.8,
    labels = freq_labs) + theme_classic(base_size = bs) + labs(color = "Frequency",
    y = "Tail-to-signal ratio", x = "Habitat") + scale_fill_discrete(guide = "none") +
    ylim(ylim)

# add contrasts
for (i in 1:2) ggf <- ggf + geom_signif(y_position = c(-9), xmin = c(i -
    diff), xmax = c(i + diff), annotation = c("*"), tip_length = 0.02,
    textsize = ts, size = 0.4, color = "gray40")

# add larger bracket
ggf <- ggf + geom_signif(y_position = c(-8), xmin = c(1), xmax = c(2),
    annotation = c("*"), tip_length = 0.02, textsize = ts, size = 0.4,
    color = "gray40")

ggfm <- plot(conditional_effects(mod, "habitat.structure:frequency.modulation"),
    plot = FALSE)[[1]] + scale_color_viridis_d(end = 0.8, labels = c("Tonal",
    "FM")) + theme_classic(base_size = bs) + labs(x = "Habitat", y = "",
    color = "") + scale_fill_discrete(guide = "none") + theme(legend.position = c(0.93,
    0.5)) + ylim(ylim)

diff <- 0.1


# add contrasts
for (i in 1:2) ggfm <- ggfm + geom_signif(y_position = c(-12), xmin = c(i -
    diff), xmax = c(i + diff), annotation = c("*"), tip_length = 0.02,
    textsize = ts, size = 0.4, color = "gray40")

ggfm <- ggfm + geom_signif(y_position = c(-9), xmin = 1, xmax = 2,
    annotation = c("*"), tip_length = 0.02, textsize = ts, size = 0.4,
    color = "gray40")

ggam <- plot(conditional_effects(mod, "habitat.structure:amplitude.modulation"),
    plot = FALSE)[[1]] + scale_color_viridis_d(end = 0.8, labels = c("Flat",
    "AM")) + theme_classic(base_size = bs) + labs(x = "Habitat", y = "Tail-to-signal ratio",
    color = "") + scale_fill_discrete(guide = "none") + theme(legend.position = c(0.93,
    0.5)) + ylim(ylim)

# add contrasts
for (i in 1:2) ggam <- ggam + geom_signif(y_position = c(-12), xmin = c(i -
    diff), xmax = c(i + diff), annotation = c("*"), tip_length = 0.02,
    textsize = ts, size = 0.4, color = "gray40")

ggam <- ggam + geom_signif(y_position = c(-9), xmin = 1, xmax = 2,
    annotation = c("*"), tip_length = 0.02, textsize = ts, size = 0.4,
    color = "gray40")


ggdur <- plot(conditional_effects(mod, "habitat.structure:duration"),
    plot = FALSE)[[1]] + scale_color_viridis_d(end = 0.8, labels = c("Short",
    "Long")) + theme_classic(base_size = bs) + labs(x = "Habitat",
    y = "", color = "") + scale_fill_discrete(guide = "none") + theme(legend.position = c(0.93,
    0.5)) + ylim(ylim)

# add contrasts
for (i in 1:2) ggdur <- ggdur + geom_signif(y_position = c(-12), xmin = c(i -
    diff), xmax = c(i + diff), annotation = c("ns"), tip_length = 0.02,
    textsize = ts, size = 0.4, color = "gray40")

ggdur <- ggdur + geom_signif(y_position = c(-9), xmin = 1, xmax = 2,
    annotation = c("ns"), tip_length = 0.02, textsize = 4, size = 0.4,
    color = "gray40")


plot_grid(ggf, ggfm, ggam, ggdur)

3.4.3 Contrasts

3.4.3.1 Frequency modulation * habitat structure

Code
cntrs <- as.data.frame(emmeans(mod, pairwise ~ frequency.modulation *
    habitat.structure)$contrasts)
names(cntrs)[3:4] <- c("l-95% CI", "u-95% CI")

coef_table <- html_format_coef_table(cntrs, highlight = TRUE, fill = viridis(10)[7])

print(coef_table)
contrast estimate l-95% CI u-95% CI
1 no_fm open - fm open 0.815 0.205 1.443
2 no_fm open - no_fm closed -3.034 -3.284 -2.788
3 no_fm open - fm closed -2.759 -3.372 -2.133
4 fm open - no_fm closed -3.849 -4.463 -3.228
5 fm open - fm closed -3.573 -3.828 -3.324
6 no_fm closed - fm closed 0.275 -0.343 0.894

3.4.3.2 Amplitude modulation * habitat structure

Code
cntrs <- as.data.frame(emmeans(mod, pairwise ~ amplitude.modulation *
    habitat.structure)$contrasts)
names(cntrs)[3:4] <- c("l-95% CI", "u-95% CI")

coef_table <- html_format_coef_table(cntrs, highlight = TRUE, fill = viridis(10)[7])

print(coef_table)
contrast estimate l-95% CI u-95% CI
1 no_am open - am open -0.821 -1.444 -0.230
2 no_am open - no_am closed -3.495 -3.749 -3.247
3 no_am open - am closed -3.933 -4.546 -3.329
4 am open - no_am closed -2.676 -3.284 -2.065
5 am open - am closed -3.111 -3.366 -2.862
6 no_am closed - am closed -0.437 -1.062 0.153

3.5 Combined results

Code
pc1 <- extended_summary(read.file = "./data/processed/regression_model_pc1.RDS",
    n.posterior = 1000, fill = viridis(10)[7], trace.palette = my.viridis,
    remove.intercepts = TRUE, highlight = TRUE, gsub.pattern = c("b_",
        "habitat.structureclosed", "amplitude.modulationam", "frequency.modulationfm",
        "durationlong"), gsub.replacement = c("", "habitat structure",
        "amplitude modulation", "frequency modulation", "duration"),
    effects = effects, print.name = FALSE, trace = FALSE, return = TRUE)

gg_pc1 <- pc1$plot + labs(x = "Degradation (PC1)")

br <- extended_summary(read.file = "./data/processed/regression_model_blur_ratio.RDS",
    n.posterior = 1000, fill = viridis(10)[7], trace.palette = my.viridis,
    remove.intercepts = TRUE, highlight = TRUE, gsub.pattern = c("b_",
        "habitat.structureclosed", "amplitude.modulationam", "frequency.modulationfm",
        "durationlong"), gsub.replacement = c("", "habitat structure",
        "amplitude modulation", "frequency modulation", "duration"),
    effects = effects, print.name = FALSE, trace = FALSE, return = TRUE)

gg_br <- br$plot + labs(x = "Blur ratio", y = "") + theme(axis.text.y = element_blank())

ea <- extended_summary(read.file = "./data/processed/regression_model_excess_attenuation.RDS",
    n.posterior = 1000, fill = viridis(10)[7], trace.palette = my.viridis,
    remove.intercepts = TRUE, highlight = TRUE, gsub.pattern = c("b_",
        "habitat.structureclosed", "amplitude.modulationam", "frequency.modulationfm",
        "durationlong"), gsub.replacement = c("", "habitat structure",
        "amplitude modulation", "frequency modulation", "duration"),
    effects = effects, print.name = FALSE, trace = FALSE, return = TRUE)

gg_ea <- ea$plot + labs(x = "Excess attenuation", y = "") + theme(axis.text.y = element_blank())

tsr <- extended_summary(read.file = "./data/processed/regression_model_tail_to_signal_ratio.RDS",
    n.posterior = 1000, fill = viridis(10)[7], trace.palette = my.viridis,
    remove.intercepts = TRUE, highlight = TRUE, gsub.pattern = c("b_",
        "habitat.structureclosed", "amplitude.modulationam", "frequency.modulationfm",
        "durationlong"), gsub.replacement = c("", "habitat structure",
        "amplitude modulation", "frequency modulation", "duration"),
    effects = effects, print.name = FALSE, trace = FALSE, return = TRUE)

gg_tsr <- tsr$plot + labs(x = "Tail-to-signal ratio", y = "") + theme(axis.text.y = element_blank())

plot_grid(gg_pc1, gg_br, gg_ea, gg_tsr, nrow = 1, rel_widths = c(6,
    2, 2, 2))

Code
estimates <- data.frame(pc1 = pc1$coef_table$Estimate, br = br$coef_table$Estimate,
    ea = ea$coef_table$Estimate, tsr = tsr$coef_table$Estimate)

names(estimates) <- c("Degradation (PC1)", "Blur ratio", "Excess attenuation",
    "Tail-to-signal ratio")

# make estimates relative to maximum estimate in data
rel_estimates <- data.frame(lapply(estimates, function(x) x/max(abs(x)) *
    0.8))

names(rel_estimates) <- c("Degradation (PC1)", "Blur ratio", "Excess attenuation",
    "Tail-to-signal ratio")

# corrplot(as.matrix(estimates), method = 'ellipse', )

signif <- data.frame(pc1 = pc1$coef_table$`l-95% CI` * pc1$coef_table$`u-95% CI` >
    0, br = br$coef_table$`l-95% CI` * br$coef_table$`u-95% CI` >
    0, ea = ea$coef_table$`l-95% CI` * ea$coef_table$`u-95% CI` >
    0, tsr = tsr$coef_table$`l-95% CI` * tsr$coef_table$`u-95% CI` >
    0)

names(signif) <- c("Degradation (PC1)", "Blur ratio", "Excess attenuation",
    "Tail-to-signal ratio")

estimates <- as.data.frame(cbind(rownames(pc1$coef_table), stack(estimates)[,
    1], stack(rel_estimates)[, 1], stack(signif)))

names(estimates) <- c("predictor", "est", "relavite.est", "sig", "response")
estimates$signif <- ifelse(estimates$sig, "p < 0.05", "p >= 0.05")


# estimates$relavite.est <- ifelse(estimates$sig,
# estimates$relavite.est, 0)

estimates$response <- factor(estimates$response, levels = rev(c("Degradation (PC1)",
    "Blur ratio", "Excess attenuation", "Tail-to-signal ratio")))

estimates$predictor <- factor(estimates$predictor, levels = c("habitat structure",
    "frequency modulation", "amplitude modulation", "frequency", "duration",
    "frequency:habitat structure", "habitat structure:duration", "habitat structure:amplitude modulation",
    "habitat structure:frequency modulation"))

ggplot(estimates, aes(x = predictor, y = response, fill = relavite.est)) +
    geom_tile() + coord_equal() + scale_fill_gradient2(low = viridis(10)[3],
    high = viridis(10)[7], guide = "none") + geom_text(aes(label = round(est,
    3), color = signif), size = 4) + scale_color_manual(values = c("black",
    "gray")) + labs(x = "", y = "", color = "P value") + theme_classic(base_size = 10) +
    theme(axis.text.x = element_text(color = "black", size = 11, angle = 30,
        vjust = 0.8, hjust = 0.8))

3.6 Combined model metadata

Code
check_rds_fits(path = "./data/processed/", html = TRUE)
priors formula iterations chains thinning warmup diverg_transitions rhats > 1.05 min_bulk_ESS min_tail_ESS seed
1 b-normal(0, 6) Intercept-student_t(3, 0.1, 2.5) sd-cauchy(0, 4) sigma-student_t(3, 0, 2.5) simo-dirichlet(1) blur.ratio ~ frequency * habitat.structure + frequency.modulation * habitat.structure + amplitude.modulation * habitat.structure + duration * habitat.structure + mo(distance_f) + (1 | location) + (1 | treatment.replicates) 20000 4 1 10000 220 (0.0055%) 0 3942.862 NA 1717938324
2 b-normal(0, 6) Intercept-student_t(3, 83.8, 34) sd-cauchy(0, 4) sigma-student_t(3, 0, 34) simo-dirichlet(1) excess.attenuation ~ frequency * habitat.structure + frequency.modulation * habitat.structure + amplitude.modulation * habitat.structure + duration * habitat.structure + mo(distance_f) + (1 | location) + (1 | treatment.replicates) 20000 4 1 10000 220 (0.0055%) 0 4940.947 8799.022 1879993832
3 b-normal(0, 6) Intercept-student_t(3, -0.3, 2.6) sd-cauchy(0, 4) sigma-student_t(3, 0, 2.6) simo-dirichlet(1) pc1 ~ frequency * habitat.structure + frequency.modulation * habitat.structure + amplitude.modulation * habitat.structure + duration * habitat.structure + mo(distance_f) + (1 | location) + (1 | treatment.replicates) 20000 4 1 10000 220 (0.0055%) 0 5197.530 3438.294 721775339
4 b-normal(0, 6) Intercept-student_t(3, -6.3, 8.7) sd-cauchy(0, 4) sigma-student_t(3, 0, 8.7) simo-dirichlet(1) tail.to.signal.ratio ~ frequency * habitat.structure + frequency.modulation * habitat.structure + amplitude.modulation * habitat.structure + duration * habitat.structure + mo(distance_f) + (1 | location) + (1 | treatment.replicates) 20000 4 1 10000 220 (0.0055%) 0 4589.723 7674.088 1096594823

 

4 Takeaways

  • Habitat structure seems to drive most of the observed degradation
  • Frequency modulation and amplitude modulation were the signal structural features that more strongly affected transmission
  • Signal features interact with habitat structure in diverse ways, sometimes increasing and sometimes decreasing degradation
  • Degradation metrics are robust to variation in signal duration

 

5 Posterior predictive checks

Code
model_list <- c(pc1 = "./data/processed/regression_model_pc1.RDS",
    blur_ratio = "./data/processed/regression_model_blur_ratio.RDS",
    excess_attenuation = "./data/processed/regression_model_excess_attenuation.RDS",
    tail_to_signal_ratio = "./data/processed/regression_model_tail_to_signal_ratio.RDS")

ndraws <- 20

for (i in seq_len(length(model_list))) {
    cat("\n\n## ", names(model_list)[i], "\n\n")

    fit <- readRDS(model_list[[i]])
    ppc_dens <- pp_check(fit, ndraws = ndraws, type = "dens_overlay_grouped",
        group = "habitat.structure")  # shows dens_overlay plot by default
    pp_mean <- pp_check(fit, type = "stat_grouped", stat = "mean",
        group = "habitat.structure", ndraws = ndraws)
    pp_scat <- pp_check(fit, type = "scatter_avg_grouped", group = "habitat.structure",
        ndraws = ndraws)
    pp_stat2 <- pp_check(fit, type = "stat_2d", ndraws = ndraws)
    pp_loo <- pp_check(fit, type = "loo_pit_qq", ndraws = ndraws)
    pp_error <- pp_check(fit, type = "error_scatter_avg_grouped",
        ndraws = ndraws, group = "habitat.structure")
    plot_l <- list(ppc_dens, pp_mean, pp_scat, pp_error, pp_stat2,
        pp_loo)

    plot_l <- lapply(plot_l, function(x) x + scale_color_viridis_d(begin = 0.1,
        end = 0.8, alpha = 0.5) + scale_fill_viridis_d(begin = 0.1,
        end = 0.8, alpha = 0.5) + theme_classic())

    print(plot_grid(plotlist = plot_l, ncol = 2))
}

5.1 pc1

5.2 blur_ratio

5.3 excess_attenuation

5.4 tail_to_signal_ratio


 

Session information

R version 4.1.0 (2021-05-18)
Platform: x86_64-pc-linux-gnu (64-bit)
Running under: Ubuntu 20.04.2 LTS

Matrix products: default
BLAS:   /usr/lib/x86_64-linux-gnu/atlas/libblas.so.3.10.3
LAPACK: /usr/lib/x86_64-linux-gnu/atlas/liblapack.so.3.10.3

locale:
 [1] LC_CTYPE=pt_BR.UTF-8       LC_NUMERIC=C              
 [3] LC_TIME=es_CR.UTF-8        LC_COLLATE=pt_BR.UTF-8    
 [5] LC_MONETARY=es_CR.UTF-8    LC_MESSAGES=pt_BR.UTF-8   
 [7] LC_PAPER=es_CR.UTF-8       LC_NAME=C                 
 [9] LC_ADDRESS=C               LC_TELEPHONE=C            
[11] LC_MEASUREMENT=es_CR.UTF-8 LC_IDENTIFICATION=C       

attached base packages:
[1] stats     graphics  grDevices utils     datasets  methods   base     

other attached packages:
 [1] ggsignif_0.6.2     emmeans_1.8.1-1    brmsish_1.0.0      cmdstanr_0.4.0    
 [5] cowplot_1.1.1      ggdist_3.2.0       brms_2.18.0        Rcpp_1.0.11       
 [9] corrplot_0.90      viridis_0.6.3      viridisLite_0.4.2  tidyr_1.1.3       
[13] ggplot2_3.4.2      warbleR_1.1.28     NatureSounds_1.0.4 seewave_2.2.1     
[17] tuneR_1.4.4        rprojroot_2.0.3    formatR_1.11       knitr_1.43        
[21] kableExtra_1.3.4   remotes_2.4.2     

loaded via a namespace (and not attached):
  [1] backports_1.4.1      systemfonts_1.0.4    plyr_1.8.7          
  [4] igraph_1.5.0.1       splines_4.1.0        crosstalk_1.2.0     
  [7] TH.data_1.1-0        rstantools_2.2.0     inline_0.3.19       
 [10] digest_0.6.33        htmltools_0.5.5      fansi_1.0.4         
 [13] magrittr_2.0.3       checkmate_2.2.0      RcppParallel_5.1.5  
 [16] matrixStats_0.62.0   xts_0.12.2           sandwich_3.0-1      
 [19] svglite_2.1.0        prettyunits_1.1.1    colorspace_2.1-0    
 [22] signal_0.7-7         rvest_1.0.3          xfun_0.39           
 [25] dplyr_1.0.10         callr_3.7.3          crayon_1.5.2        
 [28] RCurl_1.98-1.12      jsonlite_1.8.7       lme4_1.1-27.1       
 [31] ape_5.6-2            survival_3.2-11      zoo_1.8-11          
 [34] glue_1.6.2           gtable_0.3.3         webshot_0.5.4       
 [37] distributional_0.3.1 pkgbuild_1.4.0       rstan_2.21.7        
 [40] abind_1.4-5          scales_1.2.1         mvtnorm_1.1-3       
 [43] DBI_1.1.3            miniUI_0.1.1.1       dtw_1.23-1          
 [46] xtable_1.8-4         proxy_0.4-27         stats4_4.1.0        
 [49] StanHeaders_2.21.0-7 DT_0.26              htmlwidgets_1.5.4   
 [52] httr_1.4.4           threejs_0.3.3        posterior_1.3.1     
 [55] ellipsis_0.3.2       pkgconfig_2.0.3      loo_2.4.1.9000      
 [58] farver_2.1.1         utf8_1.2.3           labeling_0.4.2      
 [61] tidyselect_1.2.0     rlang_1.1.1          reshape2_1.4.4      
 [64] later_1.3.0          munsell_0.5.0        tools_4.1.0         
 [67] cli_3.6.1            generics_0.1.3       ggridges_0.5.4      
 [70] evaluate_0.21        stringr_1.5.0        fastmap_1.1.1       
 [73] yaml_2.3.7           processx_3.8.2       purrr_1.0.0         
 [76] pbapply_1.7-2        nlme_3.1-162         mime_0.12           
 [79] projpred_2.0.2       xml2_1.3.5           brio_1.1.3          
 [82] compiler_4.1.0       bayesplot_1.9.0      shinythemes_1.2.0   
 [85] rstudioapi_0.14      gamm4_0.2-6          testthat_3.1.10     
 [88] tibble_3.2.1         stringi_1.7.12       ps_1.7.5            
 [91] Brobdingnag_1.2-9    lattice_0.21-8       Matrix_1.5-4.1      
 [94] nloptr_1.2.2.2       markdown_1.3         shinyjs_2.1.0       
 [97] fftw_1.0-7           tensorA_0.36.2       vctrs_0.6.3         
[100] pillar_1.9.0         lifecycle_1.0.3      bridgesampling_1.1-2
[103] estimability_1.4.1   bitops_1.0-7         httpuv_1.6.6        
[106] R6_2.5.1             promises_1.2.0.1     gridExtra_2.3       
[109] codetools_0.2-18     boot_1.3-28          colourpicker_1.2.0  
[112] MASS_7.3-60          gtools_3.9.3         assertthat_0.2.1    
[115] rjson_0.2.21         withr_2.5.0          shinystan_2.6.0     
[118] multcomp_1.4-17      mgcv_1.8-42          parallel_4.1.0      
[121] grid_4.1.0           minqa_1.2.4          coda_0.19-4         
[124] rmarkdown_2.23       shiny_1.7.3          base64enc_0.1-3     
[127] dygraphs_1.1.1.6